Creando una sistema de Algebra Lineal

 Ariel Vallarino

 Creo clase con métodos para realizar operaciones de álgebra lineal con vectores y matrices


In [2]:
class Array:    
    #----------[ Metodo constructor y validador:
    def __init__(self, list_of_rows): 
        """ Metodo Constructor y Validador
        2. Un validador
        
        Parámetros:
        .- list_of_rows: lista de arrays que asigna al atributo DATA 
        y ademas en SHAPE define la dimensión del mismo.
        
        Ejemplo: 
            A = Array([[1,2,3],[4,5,6],[7,8,9]])
        """
        # obtener dimensiones
        self.data = list_of_rows
        nrow = len(list_of_rows)
        
        #  ___caso vector: redimensionar correctamente
        if not isinstance(list_of_rows[0], list):
            nrow = 1
            self.data = [[x] for x in list_of_rows]
        
        # ahora las columnas deben estar bien aunque sea un vector
        ncol = len(self.data[0])
        self.shape = (nrow, ncol)
        
        # validar tamano correcto de filas
        if any([len(r) != ncol for r in self.data]):
            raise Exception("Las filas deben ser del mismo tamano")
    #---------- Metodo constructor y validador ].
    
    
    #----------[ Metodo para imprimir sin función Print:
    def __repr__(self):
        """ 1.1. Un metodo para imprimir mejor...
        Metodo para imprimir un Array sin utilizar print
        """
        nrow, ncol = self.shape
        data_out  = "array(["
        data_line = "["
          
        # Obtengo Max. longitud de los valores para generar empacio entre las columnas:
        espacio = max(max([[len(str(self.data[i][j])) for i in range(nrow)] for j in range(ncol)]))                                         
        
        for i in range(nrow):
            data_line  = "[" + str(", ".join(map(lambda i: ('{:2}'.format(i)).rjust(espacio), self.data[i]))) + "]"
            if i < (nrow - 1):
                data_out += data_line + ",\n       "
            else: 
                data_out += data_line + "])"
                                
        return (data_out)
    #---------- Metodo para imprimir sin función Print ].
 

    #----------[ Metodo para imprimir con función Print: 
    def __str__(self):
        """ 1.2. Un metodo para imprimir mejor...
        Metodo para imprimir un Array utilizando print
        
        """      
        nrow, ncol = self.shape
        data_out  = "["
        data_line = "["
        
        espacio = max(max([[len(str(self.data[i][j])) for i in range(nrow)] for j in range(ncol)]))                                  
        
        for i in range(nrow):
            data_line  = "[" + str(" ".join(map(lambda i: ('{:2}'.format(i)).rjust(espacio), self.data[i]))) + "]"
           
            if i < (nrow - 1):
                data_out += data_line + "\n "
            else: 
                data_out += data_line + "]"
                                
        return (data_out)
    #---------- Metodo para imprimir con función Print ].
    
    
    #----------[ Metodo para obtener un item: 
    def __getitem__(self, idx):
        """ 3.1. Indexing and Item assignment
        Retorna un item """    
        return self.data[idx[0]][idx[1]]        
    #---------- Metodo para obtener un item ].
    
    
    #----------[ Metodo para modificar un item:
    def __setitem__(self, idx, new_value):
        """ 3.2. Indexing and Item assignment
        Modifica el valor de un item particular """               
        self.data[idx[0]][idx[1]] = new_value            
        print("Item modificado.")        
    #---------- Metodo para modificar un item ].
    
    
    #----------[ Función para crear Matriz de Ceros:
    def zeros(shape): 
        """ 4.1. Iniciar una matriz en Ceros
        Genera un Array de Ceros de la dimensión indicada
        
        Parametros
        ----------
        shape: (n,m) -> n filas x m columnas
        """        
        if isinstance(shape, tuple):
            nrow, ncol = shape
        elif isinstance(shape, int):
            nrow = ncol = shape        
        else:
            raise Exception("Parámetro no definido.")
        
        newArray = Array([[0 for i in range(ncol)] for j in range(nrow)])      
        return newArray
    #---------- Función para crear Matriz de Ceros ]. 
    

    #----------[ Función para crear Matriz Identidad:
    def eye(n,m=0,k=0):
        """ 4.2. Iniciar una matriz con Unos en la diagonal
        
        Parametros
        ----------
        n: Número de filas
        m: Número de columnas (opcional)
           Por defecto m = n
        k: Desplaza la diagonal (opcional)
           k = 0: Daigonal principal (valor por defect)
           k > 0: Diagonal superior
           k < 0: Diagonal inferior
        """
        if m == 0: m = idxy = n
        elif n > m: idxy = n
        else: idxy = m
            
        idx = 0
        idy = 0
            
        if k == 0: 
            newArray = Array([[0 if i is not j else 1 for i in range(n)] for j in range(m)])      
        else:
            newArray = Array.zeros((n,m))            
            for i in range(idxy):
                idx = idy = i            
                if k > 0:
                    idy += k 
                    if (idx < n) & (idy < m):
                        # Asigno 1 a la diagonal
                        newArray.data[idx][idy] = 1.0  
                else:                        
                    idx += abs(k)                     
                    if (idx < n) & (idy < m) & (idx >= 0):
                        # Asigno 1 a la diagonal
                        newArray.data[idx][idy] = 1.0  
                
        return newArray
    #---------- Función para crear Matriz Identidad ].
    
    
    #----------[ Metodo para generar la matriz transpuesta:
    def transpose(self):
        """ 5. Transposicion
        Calcula la Transpuesta de una matriz
        """
        nrow, ncol = self.shape
        return Array([[self.data[j][i] for j in range(nrow)] for i in range(ncol)])
    #---------- Metodo para generar la matriz transpuesta ].
    
    
    #----------[ Metodo para sumar arrays:
    def __add__(self, other):
        """ 6.1. Suma
        Suma Arrays
        """
        if isinstance(other, Array):
            if self.shape != other.shape:
                raise Exception("Las dimensiones son distintas!")
            rows, cols = self.shape
            newArray = Array([[0. for c in range(cols)] for r in range(rows)])
            for r in range(rows):
                for c in range(cols):
                    newArray.data[r][c] = self.data[r][c] + other.data[r][c]
            return newArray
        elif isinstance(2, (int, float, complex)): # en caso de que el lado derecho sea solo un numero
            rows, cols = self.shape
            newArray = Array([[0. for c in range(cols)] for r in range(rows)])
            for r in range(rows):
                for c in range(cols):
                    newArray.data[r][c] = self.data[r][c] + other
            return newArray
        else:
            return NotImplemented # es un tipo de error particular usado en estos metodos
    #---------- Metodo para sumar arrays ].


    #----------[ Metodo para sumar arrays:
    def __radd__(self, other):
        """ 6.2. Suma
        Suma Arrays
        """
        return self.__add__(other)
    #---------- Metodo para sumar arrays ].
    
    
    #----------[ Metodo para restar arrays:
    def __sub__(self, other):
        """  6.3. Resta 
        Resta Arrays
        """         
        if isinstance(other, Array):
            if self.shape != other.shape:
                raise Exception("Las dimensiones son distintas!")
            rows, cols = self.shape
            newArray = Array([[0. for c in range(cols)] for r in range(rows)])
            for r in range(rows):
                for c in range(cols):
                    newArray.data[r][c] = self.data[r][c] - other.data[r][c]
            return newArray
        elif isinstance(2, (int, float, complex)): # en caso de que el lado derecho sea solo un numero
            rows, cols = self.shape
            newArray = Array([[0. for c in range(cols)] for r in range(rows)])
            for r in range(rows):
                for c in range(cols):
                    newArray.data[r][c] = self.data[r][c] - other
            return newArray
        else:
            return NotImplemented
    #---------- Metodo para restar arrays ].
    
    
    #----------[ Metodo para sumar arrays:
    def __rsub__(self, other):
        """ 6.3. Resta 
        Resta Arrays
        """        
        return self.__sub__(other)
    #---------- Metodo para sumar arrays ].
    
    
    #----------[ Metodo para Multiplicar arrays:
    def __mul__(self, other):
        """ 7.1. Multiplicacion Matricial
        Multiplica Arrays:
        """
        if isinstance(other, Array):              
            other_t = other.transpose() # Obtengo transpuestas
            idx = 0
            idy = 0                        
            
            if "Vector" not in str(type(other)):                   
                newArray = Array([[sum(m*n for m,n in zip(j, i)) for j in other_t.data] for i in self.data])
                
                return newArray
            else: 
                list_of_calc = [0 for x in range(self.shape[0])]  # Creo lista de long. n con Ceros          
                
                for x in range(self.shape[0]):   # Itero n                 
                    # Calculo la suma de multiplicar Fila Nx por Columna My                                                               
                    list_of_calc[idx] = (sum(i*j for i,j in zip(other_t.data[idy], self.data[idx])))
                    idx += 1 
                
                return Vector(list_of_calc)
    
        elif isinstance(other, int):            
            if "Vector" not in str(type(self)):                
                nrow, ncol = self.shape
                return Array([[self.data[i][j] * other for j in range(nrow)] for i in range(ncol)])
                
            else:
                return NotImplemented
        else:
            return NotImplemented
    #---------- Metodo para Multiplicar arrays: ].
    
    
    #----------[ Metodo para sumar arrays:
    def __rmul__(self, other):
        """ 7.2. Multiplicacion Matricial
        Multiplica Arrays:
        """
        return self.__mul__(other)
    #---------- Metodo para sumar arrays ].
    
#------------------------------------------------------------------------------------------#  

    #----------[ Metodo para resolver Ly = b:
    def forward_subs(pL, pb):
        """ Funcion forward_subs que resuelva sistemas de ecuaciones de la forma: 
                Lx = y
            con L triangular inferior e y cualquier Vector o Array de una columna
        """
        pass # Valdiar dimensiones y que L sea Triangular Inferior
    
        iy = [0 for n in range(pL.shape[0])]        
        
        for i in range(pL.shape[0]):    
            iy[i] = pb.data[i][0]            
            for j in range(pL.shape[1]):                                 
                if j < i:
                    iy[i] -= (pL.data[i][j] * iy[j])
            
            iy[i] /= pL.data[i][i]
            iy[i] = round(iy[i], 3)
            
        return Vector(iy)    
    #---------- Metodo para resolver Ly = b ].
    
    
    #----------[ Metodo para resolver Ux = y:
    def backward_subs(pU, py):        
        """ Funcion backward_subs que resuelva sistemas de ecuaciones de la forma: 
                Ux = y
            con U triangular superior y y Vector o Array de una columna. 
        """
        pass # Valdiar dimensiones y que U sea Triangular Superior
    
        ix = [0 for n in range(pU.shape[0])]        
        
        for i in range(pU.shape[0], 0, -1):            
            i -= 1
            ix[i] = py.data[i][0]                        
            for j in range(pU.shape[1] - 1, i -1, -1): 
                if j > i:
                    ix[i] -= (pU.data[i][j] * ix[j])
            
            ix[i] /= pU.data[i][i]
            ix[i] = round(ix[i], 3)
            
        return Vector(ix)    
    #---------- Metodo para resolver Ux = y ].
    
    
    #----------[ Metodo para determinar P:
    def permutacion(self):     
        """ A parir de una matriz A genera matriz de permutaciones """
        
        nrow, ncol = self.shape
        # crea Matriz Identidad:                                                                                                                                                                                           
        mid = Array.eye(nrow, ncol)
       
        # Ordena Matriz identidad para que los valor mayores
        # coincidan con los elemetos de la diagonal:
        for j in range(ncol):        
            row = max(range(j, ncol), key=lambda i: abs(self.data[i][j]))
            if j != row:
                mid.data[j], mid.data[row] = mid.data[row], mid.data[j]

        return mid
    #---------- Metodo para determinar P ].
    
    
    #----------[ Metodo para descomponer A en LUP:
    def lu_decomposition(self):
        """ funcion LU que recibe un Array A y devuelva 3 arrays L, U y P tales que 
                PA = LU
            con L trangular inferior, U triangular superior y P matriz de permutacion
        """
        pass # Agregar validaciones / excepiones
    
        nrow, ncol = self.shape 
        
        # Inicializo a L como una Matriz identidad: 
        ML = Array.eye(nrow, ncol)        
        # Inicializo a U con 0's: 
        MU = Array.zeros(nrow)    

        # Genero martiz de permutacion
        MP = self.permutacion()
        # Reordeno array de acuerdo a las permutaciones:
        MPA = MP * self
        
                                                                                                                                                                                                                              
        for j in range(nrow):
            for i in range(j+1):            
                vu = sum(MU.data[k][j] * ML.data[i][k] for k in range(i))
                MU.data[i][j] = MPA.data[i][j] - vu

            for i in range(j, ncol):               
                vl = sum(MU.data[k][j] * ML.data[i][k] for k in range(j))
                ML.data[i][j] = (MPA.data[i][j] - vl) / MU.data[j][j]
                    
        return (MP, ML, MU)
    #---------- Metodo para descomponer A en LUP ].
    
    
    #----------[ Metodo para resolver sistemas de ecuaciones:
    def lu_linsolve(pA,pb):
        """ Funcion lu_linsolve que resuelve sistemas de ecuaciones 
                Ax = y 
            con A un Array y y un Vector o Array de una columna
        """
        # Descomposicion de A en P, L y U:
        MP, ML, MU = pA.lu_decomposition()
        # De acuerod a Ly = b -> Obtengo y:
        sy = Array.forward_subs(ML,pb)
        # De acuerod a Ux = y -> Obtengo x:
        sx = Array.backward_subs(MU,sy)                      
        
        pA = pA * MP
        # Retorno solución:
        return sx
    #---------- Metodo para resolver sistemas de ecuaciones ].
    
#------------------------------------------------------------------------------------------#

class Vector(Array): # declara que Vector es un tipo de Array
    """8. Vectores - Herencia 
    """
    def __init__(self, list_of_numbers):
        self.vdata = list_of_numbers
        list_of_rows = [[x] for x in list_of_numbers]
        return Array.__init__(self, list_of_rows)
    
    def __repr__(self):
        return "Vector(" + str(self.vdata) + ")"
    
    def __str__(self):
        return str(self.vdata)
    
    def __add__(self, other):
        new_arr = Array.__add__(self, other)
        return Vector([x[0] for x in new_arr.data])
    
    def __sub__(self, other):        
        new_arr = Array.__sub__(self, other)
        return Vector([x[0] for x in new_arr.data])
    
#------------------------------------------------------------------------------------------#

Uso de clase Array

Generar array:


In [3]:
A = Array([[1,2,3], [4,5,6],[7,8,9]])
A.__dict__     # el campo escondido __dict__ permite acceder a las propiedades de clase de un objeto


Out[3]:
{'data': [[1, 2, 3], [4, 5, 6], [7, 8, 9]], 'shape': (3, 3)}

In [4]:
A.data, A.shape


Out[4]:
([[1, 2, 3], [4, 5, 6], [7, 8, 9]], (3, 3))

Visualizar datos:


In [5]:
print(A)  # Muestro el contenido de A utilizando la función print


[[ 1  2  3]
 [ 4  5  6]
 [ 7  8  9]]

In [6]:
A        # Muestro el contenido de A sin utilizar la función print


Out[6]:
array([[ 1,  2,  3],
       [ 4,  5,  6],
       [ 7,  8,  9]])

Accediendo a una posición:


In [7]:
A[1,:]       # Accedo al valor de una posición


Out[7]:
[4, 5, 6]

In [8]:
A[0,1] = -3  # Modifico el valor de una posición
A


Item modificado.
Out[8]:
array([[ 1, -3,  3],
       [ 4,  5,  6],
       [ 7,  8,  9]])

Zeros


In [9]:
Z = Array.zeros((3,3))   # Genero matriz de m x n de Ceros
Z


Out[9]:
array([[ 0,  0,  0],
       [ 0,  0,  0],
       [ 0,  0,  0]])

Matriz Identidad


In [10]:
E = Array.eye(3,3)     # Genero matriz de m x n con 1s en la diagonal
E


Out[10]:
array([[ 1,  0,  0],
       [ 0,  1,  0],
       [ 0,  0,  1]])

Transpuesta


In [11]:
T = A.transpose()      # Obtengo la transpuesta de A
T
print("A:")
print(A)
print("Transpuesta de A:")
print(T)


A:
[[ 1 -3  3]
 [ 4  5  6]
 [ 7  8  9]]
Transpuesta de A:
[[ 1  4  7]
 [-3  5  8]
 [ 3  6  9]]

Operaciones entre Arrays

Suma de Matrices


In [12]:
S1 = A + T          # Sumo Matrices
S1


Out[12]:
array([[ 2,  1, 10],
       [ 1, 10, 14],
       [10, 14, 18]])

Suma de Matriz + Escalar


In [13]:
S2 = A + 10        # Sumo Matriz + Escalar
S2


Out[13]:
array([[11,  7, 13],
       [14, 15, 16],
       [17, 18, 19]])

Multiplicación de matrices


In [14]:
M1 = A * T             # Multiplico Matrices
M1


Out[14]:
array([[ 19,   7,  10],
       [  7,  77, 122],
       [ 10, 122, 194]])

Multiplicación de Matriz x Escalar


In [15]:
M2 = A * 5            # Multiplico Matriz por un Escalar
M2


Out[15]:
array([[  5, -15,  15],
       [ 20,  25,  30],
       [ 35,  40,  45]])

Vectores


In [16]:
Vector([1,2,3]).__dict__


Out[16]:
{'data': [[1], [2], [3]], 'shape': (3, 1), 'vdata': [1, 2, 3]}

In [17]:
vec = Vector([1,2,3])
vec


Out[17]:
Vector([1, 2, 3])

Suma de Vectores


In [18]:
Vector([1,2,3]) + Vector([5,-2,0])


Out[18]:
Vector([6, 0, 3])

Suma de Vector + Escalar


In [19]:
(vec + 10) - 5


Out[19]:
Vector([6, 7, 8])

Matriz & Vector


In [20]:
MV = A * vec
MV


Out[20]:
Vector([4, 32, 50])

Resolviendo sistemas de ecuaciones


In [21]:
L = Array([[1,0,0,0], [2,1,0,0],[-1.5,0.625,1,0],[-1,1.75,6.6667,1]])
b = Vector([4,-8,-4,-1])
U = Array([[2,3,2,4],[0,4,-8,-8],[0,0,3,9],[0,0,0,-49]])

Forward Subs: Ly = b


In [22]:
y = Array.forward_subs(L,b)
y


Out[22]:
Vector([4.0, -16.0, 12.0, -49.0])

Backward Subs: Ux = y


In [23]:
x = Array.backward_subs(U,y)
x


Out[23]:
Vector([-1.0, 0.0, 1.0, 1.0])

Descomposición LU


In [37]:
A1 = Array([[4,10,-4,0],[-2,4,4,-7],[-3,-2,-5,-2],[2,3,2,4]])
b1 = Vector([-8,-1,-4,4])
print("A1: ")
print(A1)
print("\nb1: ")
print(b1)


A1: 
[[ 4 10 -4  0]
 [-2  4  4 -7]
 [-3 -2 -5 -2]
 [ 2  3  2  4]]

b1: 
[-8, -1, -4, 4]

Dada una array A genero 3 arrays L, U y P tales que PA = LU


In [38]:
P1, L1, U1 = A1.lu_decomposition()
print("P: " )
print(P1)
print("\nL: ")
print(L1)
print("\nU: ")
print(U1)


P: 
[[ 1  0  0  0]
 [ 0  1  0  0]
 [ 0  0  1  0]
 [ 0  0  0  1]]

L: 
[[  1.0     0     0     0]
 [ -0.5   1.0     0     0]
 [-0.75 0.6111111111111112   1.0     0]
 [  0.5 -0.2222222222222222 -0.48192771084337355   1.0]]

U: 
[[                 4                 10                 -4                  0]
 [                 0                9.0                2.0               -7.0]
 [                 0                  0 -9.222222222222221 2.2777777777777786]
 [                 0                  0                  0 3.5421686746987957]]

Resuelve sistema Ax = b


In [39]:
solucion = Array.lu_linsolve(A1,b1)
solucion


Out[39]:
Vector([-1.0, 0.0, 1.0, 1.0])

In [ ]: